
* working folder:
	cd "....."
	
	
********************************************************************************
********************************************************************************

	*** Making all Figures:
	
	use "tax_cleanedsample_stage1", clear
	
	
	** Figure 2: sample distribution around the threshold
	
	preserve
	
	local win=20
	keep if avgsale>=25-`win' & avgsale<=25+`win'
	
	global s s1mono
	
	twoway (hist avgsale, fraction width(1) ) (scatteri 0 25 0.075 25, recast(line) lp(dash) lc(red)) , ///
		xlabel(0(5)50) xtitle(3-year average sales ($ millions)) ylabel(0(0.025)0.075) ytitle(Fraction of sample) ylabel(,grid) ///
		name(fig_firmdist, replace) scheme($s) legend(off)
	restore
	
	
	** Figure 3: McCrary test

	DCdensity avgsale if avgsale>5 & avgsale<45, breakpoint(25) generate(Xj Yj r0 fhat se_fhat) graphname(DCdensity_nolog.eps) //density test developed by McCrary (2008)
			
	restore

	***************************************************************************************************
	***************************************************************************************************
	
	*** Constructing the main RDD sample:
	
	global covariates size q_ind prof tang
	global basedepvar lev inv comstock_c totfinratio
	global othervars sale avgsale treated switching
	global othervars2 totalassetscurrent capital Lsize Lq_ind Lprof Ltang idit opincome gp ndshield_gp interest int_opincome
	
	* Collapsing the panel data to firm level before and after:
	collapse (last) naics naics2 naics4 public (mean)  $basedepvar $depvar2 $covariates $othervars $othervars2, by(newid after)

	* Keep firms that exist in the data both before and after the shock:
	by newid: gen cnt=_N
	keep if cnt==2

	sort newid after
	ds $covariates
	foreach var in $basedepvar `r(varlist)' {
		sort newid after
		by newid: gen D`var'=`var'[_N]-`var'[1]
	}
	
	foreach var in $basedepvar $covariates $othervars $othervars2 {
		sort newid after
		by newid: gen `var'_pre=`var'[1]
	}

	by newid: egen grandavgsale=mean(avgsale)
	forvalues j=1/4 {
		gen salediffP`j'=(grandavgsale - 25)^`j'
	}

	by newid: keep if _n==_N
	
	winsor2 D* $basedepvar $depvar2 $covariates $othervars $othervars2 , replace cuts(3 97)
	
	keep if ~switching
	keep if treated==0 | treated==1
	
	destring naics4, replace
	
	
	***************************************************************************************************
	***************************************************************************************************

	*** Figure 4: RDD figures for predetermined firm characteristics:
	
	global graphopts legend(pos(6) col(2) region(lstyle(none))) xtitle(Average Sales) xlabel(10 20 25 30 40, grid glpattern(dot) ) scheme(s1mono)
		
	local power 2
	local win 15
	local numbin 10
	
	preserve	
	keep if avgsale>=25-`win' & avgsale<=25+`win'
	
	* using rdplot
	rdplot q_ind_pre avgsale, c(25) p(`power') h(`win') nbins(`numbin') ci(95) ///
		graph_options($graphopts ytitle(Market-to-book) name(fig_mktbk_pre, replace) ylabel(0(2)5, grid glpattern(dot)) )
		
	rdplot prof_pre avgsale, c(25) p(`power') h(`win') nbins(`numbin') ci(95) ///
		graph_options($graphopts ytitle(Profitability) name(fig_prof_pre, replace) ylabel(-.25(.125).25, grid glpattern(dot)) )

	rdplot tang_pre avgsale, c(25) p(`power') h(`win') nbins(`numbin') ci(95) ///
		graph_options($graphopts ytitle(Tangibility) name(fig_tang_pre, replace) ylabel(0(.1).5, grid glpattern(dot)) )
	
	rdplot lev_pre avgsale, c(25) p(`power') h(`win') nbins(`numbin') ci(95) ///
		graph_options($graphopts ytitle(Debt ratio) name(fig_lev_pre, replace) ylabel(.25(.25)1, grid glpattern(dot)) )

	rdplot inv_pre avgsale, c(25) p(`power') h(`win') nbins(8) ci(95) ///
		graph_options($graphopts ytitle(Investment rate) name(fig_inv_pre, replace) ylabel(-.2(.1).3, grid glpattern(dot)) )

	restore
			
	***************************************************************************************************
	***************************************************************************************************

	*** Figure 5: RDD figures for firm outcomes:
		
	preserve

	global graphopts legend(pos(6) col(2) region(lstyle(none))) xtitle(Average Sales) xlabel(10 20 25 30 40, grid glpattern(dot) ) scheme(s1mono)
	
	local power 2
	local win 15
	local numbin 5
	
	keep if avgsale>=25-`win' & avgsale<=25+`win'

	qui reg Dlev Dsize Dq_ind Dprof Dtang i.naics4 , vce(cluster naics2)
	predict Dlev_res, residuals
	rdplot Dlev_res avgsale, c(25) p(`power') h(`win') nbins(`numbin') ci(90) ///
		graph_options(ytitle({&Delta} Debt ratio) name(fig_dlev, replace) ylabel(-.08(.04).12, grid glpattern(dot)) $graphopts )
	
	qui reg Dcomstock_c Dsize Dq_ind Dprof Dtang i.naics4 , vce(cluster naics2)
	predict Dcomstock_c_res, residuals
	rdplot Dcomstock_c_res avgsale, c(25) p(`power') h(`win') nbins(8) ci(95) ///
		graph_options($graphopts ytitle({&Delta} Common stock ratio) name(fig_dcomstock, replace) ylabel(-.3(.15).3, grid glpattern(dot)) )

	qui reg Dinv Dq_ind Dprof i.naics4 , vce(cluster naics2)
	predict Dinv_res, residuals
	rdplot Dinv_res avgsale, c(25) p(`power') h(`win') nbins(`numbin') ci(90) ///
		graph_options($graphopts ytitle({&Delta} Investment rate) name(fig_dinv, replace) ylabel(-.05(.025).025, grid glpattern(dot)) )

	restore
	
	
	
	
